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Abstract 

We investigate numerically the dynamics of unstable gravity driven three-dimensional thin liquid 
films on hydrophilic-hydrophobic patterned substrates. We explore longitudinally striped and 
checkerboard arrangements. Simulations show that for longitudinal stripes, the thin film can be 
guided preferentially on the hydrophilic stripes, while fingers develop on adjacent hydrophobic 
stripes if the width of the stripes is large enough. On checkerboard patterns, the film develops 
as a finger on hydrophobic domains, while it spreads laterally to cover the hydrophilic domains, 
providing a mechanism to tune the growth rate of the film. By means of kinematical arguments, we 
quantitatively predict the growth rate of the contact line on checkerboard arrangements, providing 
a first step towards potential techniques that control thin film growth in experimental setups. 
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FIG. 1: Schematic top view of sawtooth and finger formation in unstable driven liquid films 
on homogeneous substrates. The contact line, shown in black, destabilizes and forms sawtooth 
structures on hydrophilic substrates (a) indicated by dark grey, while fingers grow on hydrophobic 
substrates (b), shown in white. Each structure has a typical periodicity, given by the most unstable 
wavelength of the front in the linear regime, Amao.) which depends of the applied forcing, the 
inclination angle of the substrate and the wetting properties of the fluid. 

I. INTRODUCTION 



Forced liquid thin films appear in manypro cesses, such as solid coating, wetting of bio 
tissues and microfluidic flow orientation uJ-3|. For such driven films advancing on dry 

nnii 

substrates, the forcing of the film triggers the destabilization of the contact line |J-l7[, 
which gives rise to interfacial structures of different shapes, as depicted in Fig. [1], depending 
on the driving force, the wettability of the substrate, and, for gravity driven films, on the 
inclination angle of the substrate l|, ^j. On hydrophilic substrates, the instability typically 
causes a deformation of the film that has the shape of saturated sawtooth structures, such 
as the ones shown in Fig. [T] (a) , which do not grow in time, and that propagate along the 
forcing direction at the typical injection velocity U. On the contrary, if the substrate is 
hydrophobic, the contact line breaks into steadily growing fingers, as in Fig. [1] (b), which 
propagate in the direction of forcing, growing at an intrinsic rate. Both sawtooth and finger 
structures have a typical transverse periodicity, which is characterized by the most unstable 
wavelength in the linear regime, Amax- Again, this intrinsic lengthscale depends on the 



applied forcing, the inclination angle, and on the wetting properties o 
Due to its relevance in technological and biological applications. 



the substrate. 
3|, new insight on 
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techniques for handling hquid films is needed. A promising scenario is that of exploiting 
substrate heterogeneity to manipulate the motion of thin films, specially in the emerging 
field of microfluidics |8| , where a thin film geornetry is a small friction alternative to common 
channels, as a medium to conduct microfiows [9j]. Substrate heterogeneity has proven useful 
for interface handling, as in the sorting of drops running on dry substrates 10|, or the 
creation of mixing domains in microchannels 

By profiting from heterogeneous substrates, our main motivation in this work will be to 
find ways in which thin films can be handled, orienting them toward a preferred path and 
controlling the way in which contact line structures grow. We shall perform a theoretical 
study, based on full three-dimensional hydrodynamical simulations, of the effect of imposing 
a hydrophilic-hydrophobic chemical pattern on the growth of the unstable contact line, 
exploring a variety of pattern configurations to gain understanding on how to control the 
motion of the thin film, on top of the intrinsic time and lengthscales of the contact line 
instability. 

The problem of unstable contact lines in driven thin films has been actively studied during 
the last two decades. Several experimental studies have addressed the problem of films 
driven by gravity on homogeneous substrates. In general, experiments are performed either 
by forcing a constant volume of liquid or by driving the film at constant injection rate. The 
constant volume configuration has been studied by Silvi and Dussan sl , Jerret and de Bruyn 
6|, and de Bruyn j3] using dry solids. They all agree in that sawtooth shaped structures are 
obtained on hydrophilic substrates, while fingers grow when using a hydrophobic substrate. 
Veretennikov, Indeikina and Chuang 12|] observe the growth of fingers on a dry substrate, 
while for a prewet substrate sawtooth shapes are observed. For films driven at a constant 



flux, Johnson, Schutler, Miksis and Bankoff 



13| observe that sawtooth and finger structures 



can be obtained for the same fluid and substrate combination by varying the inclination 
angle. 

Theoretical efforts have focused on the early stages of contact line destabilization and 
the subsequent nonlinear dynamics on homogeneous substrates. Generally speaking, the 
framework is that of lubrication theory, in which the film dynamics are reduced to a two 
dimensional problem by profiting from the smallness of the film thickness, and which is valid 
for hydrophilic substrates. Wetting properties enter in the regularization of the spurious 
divergence of viscous dissipation at the contact line. Standard approaches to this problem 
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, by increasing the thickness of the 



16| . The main effect is contained in 



either relax the no shp boundary condition, thus allowing for slip at the solid, or include 
a precursor film just in front of the macroscopic film. Within the lubrication framework, 
results for the linear stability of the contact line show that the instability is weakened by 
decreasing the inclination angle of the substrate 14 
precursor film 15| or by increasing the slip velocity 
the three dimensional structure of the film at the proximity of the contact line, where the 
free surface is smoothed by surface tension and a capillary ridge is formed. The thickness of 
the capillary ridge decreases as the wettability of the surface increases, or as the inclination 
angle of the substrate decreases. As explained by Brenner Q], when a flat contact line is 
perturbed, spatial variations of the thickness of the ridge trigger transverse capillary flows 
from the troughs to the peaks of the perturbation. As a result, the ridge thickens at the 
peaks and the perturbation grows in time, the growth rate and band of unstable modes 
being determined by the thickness of the ridge. At long times, nonlinear structures emerge. 
Moyle, Chen and Homsy [l7| studied this regime numerically using the slip model, while 
Kondic and Diez 18| used the precursor film model. These studies indicate that whatever 
the shape of the contact line, sawtooth or finger, there exists a natural length scale for finger 
or sawtooth spacing which is roughly determined by the most unstable wavelength predicted 
by linear theory. 

On perfectly homogeneous substrates, one expects that growing fingers are equally spaced 
because of the intrinsic length scale scale imposed by the most unstable wavelength. Exper- 
iments show that this ideal configuration is not easy to obtain, as finger spacing has a high 
dispersion, of about 30% js, 7|, due to, e.g., substrate impurities. In an attempt to control 
finger spacing, chemically patterned substrates have been studied both experimentally and 

l ary 



numerically. Kataoka and Troian |2[ studied the growth of fingers driven by thermocapi 



19 



20| 



stresses on alternated low and high energy longitudinal stripes. Kondic and Diez 
performed experiments and numerical simulations of the lubrication equations of gravity 
driven films on patterned substrates composed of low and high flow resistance longitudinal 
stripes. Zhao and Marshall jsj performed numerical simulations of the lubrication equations 
considering sinusoidal variations of fluid wetting properties. Both experiments and numer- 
ical simulations show that for such longitudinal patterns, the fluid advances preferentially 
over the domains that have a smaller flow resistance. As a consequence, the spacing between 
fingers can be controlled if the width of the stripes is comparable to the most unstable wave- 
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length of the front. For thinner stripes, surface tension smooths the interface, creating wider 
structures. On the other hand, wider stripes eventually allow for the growth of unstable 
wavelengths, which generate additional fingers. 

In this paper we will study the orientation of the thin film on hydrophilic-hydrophobic 
patterns of sharply contrasted wetting properties, thus departing from the hydrophilic regime 
of previous numerical studies [s], [lol, 20|, analyzing the interplay between the typical length- 



scale of the pattern and the intrinsic lengthscale of the instability. Apart from orienting the 
film on a prescribed path, another interesting problem is how to control the growth of the 
contact line structures triggered by the instability. Both experimental and theoretical results 
show that fingering can be suppressed by either decreasing the inclination angle of the sub- 
strate or by increasing its wettability sufficiently, as in these limits the contact line saturates 
to a sawtooth shape. Decreasing the inclination angle is effective only for mostly wetting 
fluids; for non wetting fluids flngering is observed even at small inclination angles js- 7, 12 1. 



Increasing surface wettability is then an appealing strategy to improve surface coverage. A 
way to increase the wettability of the substrate is by imposing a chemical pattern on the 
solid. Such an exploitation of substrate heterogeneity to control the growth of the fllm has 
not been studied previously. This shall be a second question to be addressed in this work. 

Overall, the instability poses two challenges when trying to handle thin fllms. First, there 
is the presence of an intrinsic lengthscale, A^ax? which imposes the typical spacing between 
contact line structures. Secondly, there is an intrinsic timescale, which corresponds to the 
growth rate of the contact line. We shall therefore explore the effect of patterns where the 
typical time and lengthscales are comparable to h^max and the typical growth rate of the 
fllm, respectively. Hence, in a hydrophilic domain, the front is expected to saturate to a 
sawtooth, while for hydrophobic domains contact line growth should occur. 

The rest of this paper is organized as follows. In section [TTl we present the hydrody- 
namical model and the lattice-Boltzmann integration algorithm that we use to perform the 
numerical simulations. Section IIIll is devoted to the results of gravity driven thin fllms 
on different substrates. In section IIII Al we study the formation of sawtooth- and flnger- 
shaped fronts on homogeneous substrates. We follow by analyzing the growth of the contact 
line in different patterned substrates. Section IIV Al is devoted to longitudinal patterns of 
hydrophilic-hydrophobic stripes. In section IIVBI we consider arrangements of hydrophilic- 
hydrophobic domains in checkerboard patterns. Finally, in section |V] we discuss our results 
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and present the conclusions of this work. 



II. DIFFUSE INTERFACE MODEL FOR DRIVEN THIN FILMS 



In this section we present the diffuse interface framework on which we rely to study thin 
film dynamics on chemically heterogeneous substrates. In a recent work 2l|, we have studied 
the dynamics of thin films on homogeneous substrates in the aforementioned framework, 
which considers the hydrodynamics of two immiscible liquid phases coupled by a diffuse 
interface. Such a model allows for contact line dynamics and has been validated against 
molecular dynamics simulations 22|. Contact line dynamics emerges as a consequence of 



a diffusive mechanism at the contact line, making the use of explicit boundary conditions 



at the contact line unnecessary 
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24j . This model has been used to study the dynamics 



of the contact line between parallel plates on hydrophilic-hydrophobic domains by Wang 
et al. {25!. Lattice-Boltzmann simulations, the powerful computational fluid dynamics 
algorithm that we shall use in this work, have been used to apply this model on a variety of 
situations. Briant and Yeomans 



23| studied the motion of contact lines between shearing 



parallel plates. 



'he development of thin films in Hele-Shaw cells was studied by Ledesma- 



Aguilar et al. 2J|. Three-dimensional 



Ledesma-Aguilar et al. 



2m- In Ref. 



effects on the Saffman- Taylor finger was studied by 



21| we have performed lattice-Boltzmann simulations 



of thin films flowing on homogeneous substrates of arbitrary wetting properties. Our results 
for the linear stability of the front show that the contact line is unstable to a wider band 
of perturbation modes as the equilibrium contact angle increases. To validate our model, 
we have shown that the previous lubrication theory results are approached for sufficiently 
small driving velocities. In addition, quantitative agreement has been obtained between 
the lattice-Boltzmann thin film profiles and those obtained by Spaid and Homsy [l^ for 
the slip-velocity lubrication model in the same small-velocity regime. Long time dynamics 
show that a sawtooth saturates when the fluid wets the substrate, while fingers grow when 
forcing non-wetting fluids. In agreement with previous results 



13 



18| , we have observed that 



sawtooth shaped fronts are obtained by decreasing the inclination angle of the substrate. 
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Hydrophobic 



H 



FIG. 2: Schematic representation of the system geometry. A thin film of thickness, h^-, width 
W , and initial length, Lq, is forced by gravity in the x direction on an heterogeneous substrate 
composed of hydrophilic {9e = 0°) and hydrophobic {9e = 90°) domains. The film is surrounded 
by a fluid of much smaller viscosity. The total system occupies a box of size L x W x H . 

A. Governing equations 

In this section we present the diffuse interface model that governs the dynamics of a 
gravity driven thin film in contact with a solid substrate. Figure [2] shows a schematic rep- 
resentation of the system, which is composed of two immiscible liquid phases of different 
viscosities in contact with a solid wall. One phase corresponds to the driven thin film, while 
the other plays the role of a surrounding fluid of much smaller viscosity. Instead of con- 
sidering that these phases are separated by a sharp interface, we introduce a concentration 
variable, or order parameter, 0(r), which has constant values in the bulk of each phase, 
= 1 for the thin film, and = — 1 for the surrounding fluid. The concentration varies 
smoothly between these values across the interface. The dynamics of both phases follows 
by imposing that mass, momentum, and concentration are conserved in time. 

For an incompressible fluid, mass conservation gives the continuity equation 



where v(r) is the velocity of the fluid. Momentum conservation leads to the Navier-Stokes 
equations 



The left hand side of equation ([21) contains the usual acceleration and inertial terms of 
the velocity field. The right hand side is composed of the usual pressure gradient term, 
where P is the pressure field, the body force term, where p is the fluid density and g is the 
acceleration due to gravity, and the viscous friction term of a Newtonian fluid, where r){4>) 



V-v = 0, 



(1) 
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is the viscosity of the fluid, which depends on the order parameter. The additional term 
in equation ([2]) is the chemical force per unit volume of fluid, associated with the gradient 
of the chemical potential /i(0). This force is important only where the concentration proflle 
varies signiflcantly, i.e., at the diffuse interface. In the sharp interface limit this interfacial 



forcing gives rise to the Young-Laplace condition (see Ref. (27j for details). 

The dynamics of the order parameter is given by the conservation equation 9t0 + v- V0 = 
— V-j, where the current, j, is given by j = — MVyU. Hence, (p obeys the convection diffusion 
equation 

^ + V • V0 = MVV, (3) 

where M is a mobility. 

In equilibrium, the state of the system is determined by the free energy functional 

p] = Jdv p) + ^(V0)2) . (4) 

The flrst term in equation (|4]) is a volume contribution, V{(f), p) = AcfP' /2 + 50^/4 + p/3 Inp, 
which allows for two phase coexistence via the dependent terms, and contains an ideal gas 
contribution given by the p dependent term. The second contribution to the free energy is 
the squared gradient term, which accounts for the energy cost of spatial variations of the 
order parameter by a factor k. From this energy functional one can work out expressions for 
the chemical potential, p, the total pressure tensor 7-*, the surface tension of the fluid-fluid 
interface, 7, the interface width, and the equilibrium values of the order parameter in the 
bulk of the phases (p* = ^(peq- The chemical potential is given by, 

p = ^ = A(l) + B(t)^ - KV^(f). 

0(p 

In this expression, p corresponds to the energy cost of varying the concentration (p, hence, it 
should be understood as the exchange chemical potential of the system. The total pressure 
tensor, T*, is composed of a density dependent contribution arising from the dependence of 
V on p, and a "chemical" contribution, which depends on the order parameter, 

V = P{p)S + + - ^ + i| V0|2)) S 

where P = ^p and 5 is the diagonal matrix. Thus, by including the p dependent term in 
equation (jl]) one recovers the behavior of an ideal gas in the case of a single homogeneous fluid 
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in equilibrium, characterized by vanishing order parameter gradients. Taking the divergence 
of the pressure tensor gives the force density acting on the fluid 

- V ■ P = - VP - 0V/i, 

which clarifies the origin of the chemical term in equation ([2]). In this paper we will consis- 
tently use B = —A, from which it follows that the remaining equilibrium properties of the 
system are given by 7 = {—8kA/9)^ , ^ = {—k/2A)^ and 0* = ±1. 

B. System geometry and boundary conditions 

We choose a system geometry to mimic the constant flux configuration of gravity driven 
thin films, as shown in figure [2J We choose a rectangular domain of linear dimensions 
L X W X H. Initially, the thin film occupies the volume Vq = Lq x W x he- The thin film 
is forced in the x direction while the surrounding fluid is left to evolve passively. This is 
achieved by setting the gravitational force as pg = pgx\{4> + 1)^- Similarly, the viscosity is 
varied using the mixing rule 77(0) = (77) (1 + 6ri(f)), where {rj) is the mean viscosity and 6rj 
is the viscosity contrast. The viscosity of the thin film is then rj^i = {r]){l + 6r]), while the 
viscosity of the surrounding fluid is r]_i = (77) (1 — 6r]). 

To complete the model, boundary conditions for v and have to be provided. We 
choose the system geometry so that the solid substrate is parallel to the x — y plane and 
located at 2 = as shown in figure [2l The impenetrability of the solid is ensured by fixing 
Vz{x,y,z = 0) = 0, while stick boundary conditions are imposed at the wall, Vx{x,y,z = 
0) = Vy{x,y,z = 0) = 0. We fix periodic boundary conditions in the y direction, i.e., 
v(x, y = W, z) = v(x, y = 0, z), y = W, z) = y = 0, z). In the x direction a constant 
mass inflow is kept by fixing dxv{x = {0, L}, y, z) = and dx(p{x = {0, L}, y, z) = 0. At the 
upper boundary, a shear free boundary condition is imposed by fixing dz\'{x, y,z = H) = 
0, while a vanishing concentration gradient normal to the boundary is fixed by imposing 
d;,(j){x,y,z = H) = 0. 

In equilibrium, a drop of fluid sits on the solid describing the shape of a spherical cap 
that intersects the solid boundary with an equilibrium contact angle 6'e(for vanishing 6^ the 
spherical cap is replaced by a film that covers the whole substrate). The equilibrium contact 
angle is determined by Young's Law — 75L+1 — 7 cos 6e + 75l_i = 0, where •jsLi is the surface 
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tension between the solid and fluid i. Hence, to include wetting effects in the model, it is 
necessary to account for the surface energies of the solid-fluid boundaries. Within the diffuse 
interface formulation, we use a Cahn surface free energy Fs[(f)s] = /5(0s'(i"))d5', which is 
the integral along the solid-fluid surface, S, of the free energy per unit area, fs = C(f)s, that 
depends on the local value of the order parameter at the boundary, (ps- In this expression, 
the parameter C can be varied to obtain a prescribed equilibrium contact angle. Minimizing 
the overall free energy gives the boundary condition at the solid 

1 dfs 



— {x,y,z = 0) , 
oz K a(ps 



(5) 



from which the equilibrium contact angle is related to the model parameters by the expres- 
sion |28| 



cos 9p 



The desired chemical pattern is obtained by imposing equation to fix the gradient of the 
order parameter at the wall. 

The diffuse interface model naturally allows for local slip at the contact line region. 



As noted by Briant and Yeomans 



22] . even though a stick boundary 



and Qian et al. 

condition is imposed to the velocity field, the diffusive term in equation ([3]) allows for motion 
of the contact line. Such diffusive flux has been studied, e.g., by Denniston and Robbins 
29[ in miscible displacements of binary fluids using molecular dynamics. In our model, the 
diffusion mechanism induces relaxation of the stick boundary condition in a length scale Id- 
Close to the contact line, at scales comparable to Z^, one observes slip, while for scales much 
larger than this lengthscale stick is recovered. Detailed studies of this length scale have been 



performed by Briant and Yeomans 



23| and Qian et al. 



22|. 



C. Dimensionless numbers and units 



In order to relate our simulation units with experimental ones, let us briefly review the 
limit of the model in the lubrication regime. 

As explained before, in the sharp interface limit the usual continuity and Navier-Stokes 
equations are recovered from equations ([T]) and ([2]), while the chemical potential term in 
equation ([2]) gives the Young-Laplace condition at the interface, 

AP = -7/C, 

10 



where AP is the pressure jump at the interface and K, is the interface curvature. 

Due to the smallness of the film thickness, he, in the lubrication limit only the in-plane 
components of the velocity field, and Vy, are considered, and velocity gradients are as- 
sumed to occur primarily in the perpendicular direction, z. Additional assumptions are that 
the flow takes place at small velocities and that it is stationary. Hence, the left hand side 
of equation (|2]) can be neglected. Following these assumptions, it is possible to average the 
flow field in the z direction to give 

(v) = --^(7VV^/i-pg), (6) 

where (v)(x, y) is the two dimensional average velocity field, and h{x, y) is the local thickness 
of the film. In this equation, the major contribution to the local pressure is capillary, 
P = -7/C - -7V2/1. 

In this limit, it follows that the coordinates, velocity, and time in equation can be 
rescaled as 

(v)- = M 1, (7) 

using units Xc = hc{3Ca)~^ , U = hlpg,j. /{3ri j^i) and tc = xJJ'^. In these expressions, the 
capillary number, Ca, measures the ratio between viscous and capillary forces, and is defined 
as Ca = ?7+if//7. 



D. Lattice-Boltzmann algorithm 

To integrate equations ([T]), ([2]) and (E]) numerically, we use LUDWIG, a lattice-Boltzmann 
parallel implementation for binary fluids {28 1. Space is represented by a lattice of nodes that 
are connected by links. The set of links determines a set of velocity vectors {cj}. Here we 
choose the cubic lattice D3Q19 model for the velocity set {cj}, which consists of nineteen 
velocity vectors (eighteen pointing towards nearest, next nearest, and next to next nearest 
neighbors) and one accounting for rest particles. The lattice-Boltzmann method consists of 
the integration of two linearized discrete Boltzmann equations 

Mr + c„ t + 1) - /,(r, t) = 1 (/,(r, t) - /^(r, t) + F/^ , (8) 



11 



and 



gi{v + Q, t + 1) - ^?,(r, t) = - {gi{v, t) - g^%r, t)) , (9) 



^3 



where fi and (^j are velocity distribution functions of the set of allowed velocities {cj}. In 
these expressions, lattice units for length, time, and mass are fixed to Ax = 1, At = 1, and 
Am = 1, respectively. 

The dynamics of the distribution functions consists of two steps. First, the distribution 
functions undergo a collision step, which corresponds to the right hand side of equations ([8]) 
and iQ, where they are relaxed to equilibrium distribution functions, denoted by f^"^ and 
g^*^. The relaxation timescale of this process is given by the parameters r/ and Tg. We fix 
Tg = 1, while Tf is related to the fluid viscosity through r/ = {Grj + l)/2. The extra term 
F/ corresponds to the external gravitational field 30|. After the collision step, the fi and 
gi are propagated to the neighboring nodes of the lattice according to the left hand side of 
equations ([8]) and ([9]). 

The mapping between the lattice-Boltzmann algorithm and the diffuse interface model 
presented above follows from the definition of the moments of the distribution functions. 
The fi are related to the fluid density and momentum by p = J2i fi pv = J2i fi^ii 
while the gi are related to the order parameter by = gi. The relaxation process ensure 
the conservation laws Y.ifi = Y.ifT^ Y^ifi^i = T^ifi^^i^ and ^.gi = Y.i9?- In the 
simulation runs we ensure that the fluid velocities are always smaller than the speed of 
sound, Cs = 1/^3; hence the flows can be considered incompressible. To complete the 
model one needs expressions for the equilibrium distribution functions and for the forcing 
term in equations ([H]) and These expressions are obtained by expanding f^'^ and gl'^ , as 
well as -F/, in powers of the velocity, and then by relating their higher order moments to the 



equilibrium properties of the diffuse interface model 28|. For our implementation the /f'', 
gl'^ and F/ read, 

f-'' = puy {^l + 3v ■ Ci + ^vv : CiC, - ^v"^ + : CjCj^ , 

/ 9 3 

= pcu^ ( + 3v • c, + -vv : CiC, - -x? + : c^c, 



and 



4a;,(l- — )(f-c,(l + vc,)-v-f) 
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Here, u stands for the three possible magnitudes of the Cj set. Coefficient values are ujq = 2/9, 
LOi = 1/9 and = 1/72; = 9/2 - 7/2TrP, A{ = A^^ = l/pTiV and = 9/{2p)V - 
3dTiV; Af, = 9/2 - 21/2M/i, A{ = A^^ = 3Mii/p and = 9/(2p)M/i(l - S), where 1 is 
the unit matrix. 

It can be shown that equations ([H]) and converge to the hydrodynamic equations ([T]), 
([2]), and ([3]) by performing a Chapman- Enskog expansion. This procedure has been carried 



out elsewhere 



3l|. 



Boundary conditions in the lattice-Boltzmann scheme are implemented by fixing the 
distribution functions at the boundaries. To implement stick and no fiow conditions at 
the solid, bounce back rules are implemented at the solid nodes js^. Periodic boundary 
conditions are obtained as usual, e.g., for the /j at y = W, fi{x, y = W, z) = fi{x, y = 0, z). 
Shear free and vanishing concentration gradient conditions are obtained in a similar way. 
The Cahn wetting boundary condition given by equation ([5]), is obtained by extrapolating 
the order parameter to the wall nodes in order to fix a prescribed gradient in the normal 
direction to the wall (281]. 

Throughout this paper, the diffuse interface model parameters are fixed to the following 
values (in lattice units): p = 1, {?]) = OA, 5r] = 0.9, 7 = 0.0023, ^ = 0.57, and M = 16. 
These values generally ensure numerical stability in the simulations. The interface thickness, 
£ = 0.57, has been tested previously and proved not to give rise to important lattice artifacts 
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33| . We will consider two equilibrium contact angles, 6e = 0°, and 6^ = 90°, which 



correspond to C = 0.00167 and C = respectively. A suitable value for the film thickness, 
which is small enough to keep reasonable computational costs and large enough to ensure 
the relaxation to the volume values of the order parameter needed for proper hydrodynamic 
resolution, is he = 8. The mean velocity of the film is fixed to U = 0.005 by choosing the 
gravitational acceleration as = 3?7+if//(p/i^). Using these values, the capillary number is 
Ca = 0.41. 



III. HOMOGENEOUS SUBSTRATES 

The intrinsic time and lengthscales in forced thin films are set by the film growth rate, 
L{9e), and by the most unstable wavelength of the linked to the contact line instability, 
Amax(6'e), which depend on the wetting properties of the substrate, determined by the equi- 
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FIG. 3: The growth rate, w, of a transverse perturbation to the contact line as a function of the 
perturbation mode, k, on homogeneous substrates with different wetting properties obtained from 
lattice-Boltzmann simulations. The growth rate has a maximum that depends on the equihbrium 
contact angle, 9^. The maximum growth rate is larger if the substrate is hydrophobic. The 
corresponding most unstable wavelength, Amax(^e) = 27r/A:max; is smaller the more hydrophobic is 
the substrate. The error bars correspond to the standard deviation of the fit performed to estimate 
uj from the time evolution of the contact line amplitude. 

librium contact angle, Of,. In order to explore the effect of substrate heterogeneity on the 
dynamics of the film, we must first explore the dynamics on homogeneous substrates, to 
characterize L and h^ax- In this section we describe the dynamics of thin films on two 
separate homogeneous substrates with equilibrium contact angles 9^. = 0° and 6^ = 90°, 
respectively. Our aim is to characterize the linear instability by measuring the dispersion 
relation on each substrate, and to characterize the base structures that emerge at long times. 



A. Linear Stability 

A numerical linear stability analysis of the front can be performed by following the evo- 
lution of an initially flat contact line that is perturbed transversely with a single mode k/kc, 
where kc = 27r/xc- At short times, the perturbation either grows or decays exponentially 
with a growth rate, u{k). We follow the evolution of single mode perturbations of small 
initial amplitude Aq/W ~ 10~^, and focus on timescales where the amplitude shows an ex- 
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(b) 

FIG. 4: Evolution of unstable fronts in substrates of different wettabilities that originate from a 
slightly perturbed front at t/tc = 0- (a) Saturated sawtooth: Of, = 0° and W/xc = 13. (b) Growing 
finger: 9^ = 90° and W/xc = 8. The width of each system is chosen to match the most unstable 
wavelength as measured from figure [3l depending on the wetting properties of the substrate. The 
plot shows the contact line position at time intervals of At/tc = 13.4, starting at t/tc = 0. For the 
last time, t/tc = 134, the full three dimensional interface is shown. 

ponential dependence on time to extract the growth rate of the perturbation. The resulting 
dispersion relations for each front are shown in figure [31 The band of modes over which the 
front is unstable is much narrower for 6c = 0° than for 6c = 90°. This is due to the thin 
film structure near the contact line, which tends to accumulate more mass when the fluid is 
in contact with hydrophobic substrates than with hydrophilic ones. As a result, fronts on 
hydrophilic substrates are stable to longer wavelengths and have smaller growth rates. In a 



recent work 



2l| , we have validated our numerical method by comparing the linear stability 



results of the diffuse interface model to those of lubrication theory. A direct comparison 
is only possible for small dynamic angles and small capillary number, Ca. Given that our 
model allows for slip at the contact line, it is natural to compare with the results obtained 



by Spaid and Homsy 



16| , for which the contact line is allowed to slip by imposing a Navier 



boundary condition at the solid. As explained in the aforementioned paper, our results 
approach the lubrication theory results for small Ca. 

B. Sawtooth and Finger formation 



Figure |3] shows that the most unstable wavelengths are A^axi^^e = 0°)/xc ~ 14 and 
Amax(^e = 90°)/xc ~ 8. Thcsc wavcleugths set the typical length scale of the contact line 
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shape at long times. Taking into account this periodicity, in the following we will reduce 
computational costs by considering system sizes whose width is of the order of Amax(^e)- 
Figure H] shows the long time evolution of the most unstable wavelengths for 9^, = 0° and 
6e = 90°. We observe that after a transient, the thin film saturates to a sawtooth shape 
on the hydrophilic substrate, while on the hydrophobic substrate a finger grows at a steady 
rate. Given that the growth rate of the finger is constant, the substrate is not covered 
evenly, and a dry region is left to be covered by the trailing edge. This contrasts with the 
dynamics in the hydrophilic substrate, on which the sawtooth saturates to a finite amplitude 
and propagates as a whole with a velocity dictated by the injection rate. Intuitively, one 
can understand the reason why a finger grows and a sawtooth saturates by examining the 
structure of the thin film. For unstable fronts, the amplitude of a perturbation grows due 
to an uncompensated distribution of the thickness of the ridge. As the perturbation grows, 
the contact line becomes increasingly curved, thus increasing the strength of the restoring 
capillary force. If the thickness of the capillary ridge is small, as on hydrophilic substrates, 
the contact line can reach a sufficiently large curvature, and the distribution of the ridge 
thickness is balanced by surface tension. Thus, the contact line saturates to a sawtooth 
shape. In contrast, if the thickness of the ridge is large, as on hydrophobic substrates, 
the curvature of the contact line is always insufficient to balance the driving force, and 
fingers grow at a rate dictated by the distribution of the ridge thickness along the front. 
We have demonstrated this effect previously j2l|, where we have shown that the growth 
rate of any sawtooth or finger can be collapsed to a universal curve, as a function of the 
difference of the thickness of the ridge between the leading and trailing edges of the contact 
line, independently of the capillary number, or the wetting properties or inclination angle of 
the substrate. Static and dynamic fluid rivulets are known to be unstable to wave lengths 
which decrease with substrate hydrophobicity, gravity and fluid flow 34|-l36l|. We have never 



observed the destabilization of a finger for the regimes we have explored, but we cannot rule 
it out; such a possibility requires a more systematic analysis. 

To characterize the flow on homogeneous substrates, we compute the local mean velocity 
by averaging the velocity field, v, in the z direction: (v) = 1/h jj'dzv. Figure [5] shows plots 
of {v)h, which is a measure of the local mass flow, for both hydrophilic and hydrophobic 
substrates. The top panels in figure [5] show the x component of the flow. For the sawtooth, 
we observe only a slight increase of the flow near its tip. This is a signature of the balance 
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FIG. 5: Grey scale maps of the average flow within the film for a steady sawtooth [9e = 0°) 
and a finger [9^ = 90°) at t/tc = 134, starting from a slightly perturbed front at t/tc = 0. The 
width of the system is fixed to W/xc = 13 for the hydrophilic substrate, and W/xc = 8 for the 
hydrophobic substrate. Top panels: the flow in the x direction, {v)xh/{Uhc), is homogeneous far 
from the contact line, and increases sharply on the body of the finger. The sawtooth shows a small 
increase, due to the variation of the film thickness. Bottom panels: the lateral flow, {v)yh/{Uhc), 
shows clearly that both the flnger and the sawtooth are fed by transverse flows. For the flnger, 
streamlines are computed close to the trailing edge and show that particles flow to the body of the 
finger as they approach this region. 

between the driving force caused by the variation of mass along the contact line, and the 
restoring surface tension. For a finger, this balance cannot be sustained, and a sharp flow 
increase in the body of the finger is observed. As shown in the bottom panels in figure|5l there 
exist two stationary regions that feed both the sawtooth and the finger. For the sawtooth, 
these regions extend along the whole contact line structure sustaining its stationary shape. 
For a finger, they are located at the trailing edge of the contact line and fix a constant flow 
of mass that causes the finger to grow steadily (see, e.g., the streamlines depicted at the 
lowermost panel of the figure). 

IV. HETEROGENEOUS SUBSTRATES 

In the previous section we presented results for sawtooth and finger growth on homo- 
geneous substrates. In both hydrophilic and hydrophobic substrates, the lengthscale of 
instability is characterized by the most unstable wavelength of the linear regime, Amax{de), 
which sets the typical spacing of both sawtooth and finger structures. 
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FIG. 6: Contact line growth on a pattern composed of hydrophilic-hydrophobic stripes oriented 
in the direction of flow. The dark and light stripes correspond to the hydrophilic (grey) and 
hydrophobic (white) substrates, respectively. Their width, w, equals the most unstable wavelength 
observed for the hydrophobic substrate, Amax{de = 90°)/xc ~ 8. The first five panels show the 
contact line at time intervals of At/tc = 33.5. The first panel shows the perturbed contact line 
at t/tc = 0. The last three panels show the effect of the longitudinal pattern by comparing the 
contact line on the patterned substrate to the heterogeneous and homogeneous substrates with the 
same initial condition at t/tc = 134. Wetting properties are fixed to 6e = 0° at the grey zones, 
while at the white zones 9e = 90°. 

In this section we explore the dynamics of the thin film on different hydrophilic- 
hydrophobic heterogeneous substrates. Our aim is to examine the effect of a hydrophilic 
pattern imposed on a hydrophobic substrate on the dynamics of the film. Furthermore, 
taking into account the periodicity of the front, we will focus on patterns of characteristic 
lengthscales which are comparable to Amax{(^e = 90°). 



A. Longitudinal stripes 

In many situations, such as in microfiuidic flow guiding and filamenting {2, 9|, it is 
desirable to orient the growth of the contact line on a prescribed direction. For substrates 
with spatially varying wetting properties, the flow is expected to orient toward the more 
hydrophilic regions due to their lower flow resistance. Meanwhile, the front is expected to 
destabilize on the hydrophobic stripes if their width is of the order of the most unstable 
wavelength, Amax{(^e = 90°). Hence, our main interest will be to examine the effect of the 



18 



longitudinal pattern for different stripe widths. In this section we analyze the growth of the 
contact line on longitudinal patterns of hydrophilic-hydrophobic properties. 

We first impose a pattern of alternating hydrophilic-hydrophobic stripes of equal width w, 
oriented in the longitudinal direction. Given that the periodicity of the front is determined 
by Amax/xc ~ 8, we consider the effect of increasing the stripe width from A^ax, choosing 
the widths w/xc = 8, w/xc = 11, and w/xc = 13. Due to the periodicity of the front, 
larger widths corresponding to w > 2Amax are not explored. As in the homogeneous case, 
we explore the evolution of only one mode, for which the corresponding wavelength is fixed 
as A = w. 

Figure shows the evolution of the front for stripes of width w/Xc = 8. Given that for 
6e = 90°, Amax/xc ~ 8, the front destabilizes as it moves over the hydrophobic stripe and 
grows steadily. In contrast, when in contact with the hydrophilic stripe, the front bends and 
saturates to a sawtooth shape. The last three plots in figure [6] show the contact line for the 
heterogeneous and homogeneous substrates at the same time. By comparing these plots we 
notice that for the same injection rate the length of the contact line on the hydrophobic stripe 
is smaller than the length of the finger that would have grown on a homogeneous substrate. 
In addition, its tip is located at a less advanced position. In contrast, the sawtooth advances 
more rapidly on the heterogeneous substrate than on the homogeneous one. 

To examine this partial inhibition of the growth on the hydrophobic stripe caused by the 
longitudinal pattern, we compute the y component of the flow, which we show in figure [Tl^a) 
for w/xc = 8. The grey scale in this plot matches that of figure |5l which corresponds to 
a homogeneous substrate. The flow field exhibits essentially the same regions of transverse 
flow as in the homogeneous case. These regions are stationary and only translate forward 
as time proceeds. Therefore, the inflow to each stripe is constant in time. Nonetheless, by 
comparing the intensity of the regions between figure |5] and figure [Tl^a) one can appreciate 
that the inflow to the hydrophobic stripe is smaller for the heterogeneous substrate, while 
there is a clear increase in the flow to the hydrophilic stripe, which makes the sawtooth 
propagate with a larger velocity compared to the homogeneous case. 

This effect is what causes the contact line on the hydrophobic stripe to have a smaller 
growth rate compared to the homogeneous situation. Nonetheless, given that transverse 
flows originate in the vicinity of the boundary between adjacent stripes, the effect should 
be observable only when the stripe width, w, is small enough. Figures [7|^b), [7(c), and[71|^d) 
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show flow maps in the y direction corresponding to w/xc = 8, w/xc = 11 and w/xc = 13 at 
the early growth stage. On the hydrophobic stripes, contact line growth is inhibited more 
effectively for narrow stripes, as can be seen from the intensity of the mass supplying regions 
depicted in the figure. For w/xc = 13, the effect is almost suppressed. The inflow to the 
hydrophobic stripe thus becomes almost equal to that of the homogeneous case. For much 
wider stripes, where the front is expected to break into several fingers, only those growing 
next to the hydrophilic stripes are expected to decrease their growth rate, while those located 
far from these boundaries should grow essentially as if they were in a homogeneous substrate. 
Anyhow, the film propagates faster on the hydrophilic stripes, so partial orientation of the 
film is achieved. 

Having characterized the dynamics on stripes of equal widths, it is straightforward to 
explore patterns where the stripes have disparate widths. In figures [Ht^a) and Et^b) we show 
the evolution of the contact line on a very thin hydrophilic (resp. hydrophobic) stripe that is 
placed next to a wide hydrophobic (resp. hydrophilic) stripe. In both cases, the wide stripe 
has a width equal to Amax, while the thin stripe is small enough to avoid the growth of the 
contact line. Figure [8]^a) shows that, in contrast with the case of wide neighboring stripes, in 
this limit the substrate structure perturbs weakly the contact line growth, by increasing the 
velocity of the trailing edge because of the low resistance offered by the hydrophilic stripe. 
In contrast, figure [Ht^b) shows that the profile has almost converged to the sawtooth expected 
on a homogeneous hydrophilic substrate, but the trailing edge in the thin hydrophobic stripe 
decouples from the trailing edge of the sawtooth, and advances with a much smaller velocity. 
As a result, the film is oriented selectively on the hydrophilic channel, given that the width 
of the hydrophobic stripe remains much smaller than A^ax- This limit is precisely the one 



reported experimentally by Kataoka and Troian 2| , and by Kondic and Diez 19|, |20| and 



Zhao and Marsall [3| . A particularly interesting feature of this regime is that one can regard 
the growth of the contact line as a fingering process, which emerges as a consequence of 
the longitudinal pattern, and not because of the linear instability of the contact line. As 
such, it is possible to control the width of the emerging fingers by tuning the width of the 
hydrophilic stripes. For large enough hydrophilic stripe widths, the contact line is expected 
to deform, giving rise to sawtooth structures. However, this is expected to modified the 
morphology of the front only, and not to alter the fingering imposed by the longitudinal 
pattern. 
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FIG. 7: Grey scale maps of the y component of the local mass flow for a sawtooth and a neighboring 
finger in a pattern composed of longitudinal stripes, (a) Steady configuration at t/tc = 134 for 
stripes of width w/xc = 8. Early finger growth at t/tc = 33.5 for (b) w/xc = 13, (c) w/xc = 11 
and (d) w/xc = 8 for an initial condition consisting of a slightly perturbed contact line with a 
transverse perturbation mode k = Itx jw. The effect of the striped pattern is to decrease the finger 
growth rate on the hydrophobic stripe as the width of the stripe decreases. 
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FIG. 8: Contact line growth in wide/thin longitudinal stripe arrangements of hydrophilic (grey) 
or hydrophobic (white) properties, (a) For a wide hydrophobic stripe in contact with a thin 
hydrophilic stripe the contact line destabilizes on the hydrophobic domain, while the trailing edge 
advances on the hydrophilic domain, (b) On a wide hydrophilic stripe, the contact line forms 
a sawtooth, while the trailing edge of the contact line advances much more slowly on the thin 
hydrophobic domain. Wetting properties are fixed to = 0° at the grey zones, while at the white 
zones Qe = 90°. 

B. Checkerboard pattern 

Controlled tuning of the growth rate of the film is an appealing technique, that could 
be used, for instance, to control the filamenting of the film in chemical networks, where 
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FIG. 9: Geometry of the asymmetric checkerboard pattern. The pattern is composed of alternated 
hydrophihc (grey) and hydrophobic (white) domains. In general, the width of the upper and 
lower domains is different, and the width of the lower domains is large enough to allow for the 
growth of unstable wavelengths that give rise to contact line growth on the hydrophobic domains 
and saturation to a sawtooth shape on the hydrophilic domains. The lengths of hydrophilic and 
hydrophobic domains, and can be chosen to fix a given fraction of hydrophilic covering. 
Wetting properties are fixed to 6e = 0° at the grey zones, while at the white zones 9e = 90°. 
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FIG. 10: Contact line growth on asymmetric checkerboard patterns. Patterns are fixed by choosing 
^h/^c = 1-3 and (a) fcB = 0.25, and (b) fcB = 0.62. The first plot corresponds to a transverse 
perturbation of the contact line at t/tc = 0. Subsequent plots are taken at time intervals of 
At/tc = 67. At the bottom of (a) and (b) we show the finger that grows in a homogeneous 
hydrophobic substrate for the same time as the last contact line snapshot of the checkerboard 
pattern. By comparing (a) and (b), it is observed that the length of the growing contact line, L, 
is smaller if the fraction of hydrophilic domains is small. 
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FIG. 11: Contact line growth on a symmetric checkerboard pattern. The length and width of the 
domains is fixed to Amax{(^e = 90°). The first panel shows the initial perturbation to the front at 
t/tc = with a perturbation wavelength that equals Amax{(^e = 90°). The next panels show the 
contact line at time intervals of At/tc = 26.8. For both upper and lower domains, the contact line 
grows as a finger when in contact with a hydrophobic domain, and spreads when in contact with 
a hydrophilic domain. Overall, the net growth rate of the contact line achieves a constant value. 

different fluid filaments move at different velocities depending on the growth rate on a 
specific chemical stripe. The results presented in the previous section indicate that a thin 
film can be oriented along an imposed path by using longitudinal striped patterns, thus 
favoring filamenting. In addition, such patterns are useful to control the spacing between 
emerging filaments. In this section we propose a way to control the growth rate of the 
film by studying the effect of a pattern consisting of hydrophilic and hydrophobic domains 
arranged in an asymmetric checkerboard configuration, as displayed in figure O Such a 
pattern converges to the longitudinal pattern the limit where the length of the hydrophilic 
or hydrophobic domains, Ih and l^, is large enough, thus giving us the possibility of direct 
comparison to results presented in the last section. 

The checkerboard pattern introduces a transverse spreading mechanism, arising every 
time that the film comes into contact with a hydrophilic domain, and that appears as a 
possible means to control the motion of the contact line. As before, we are interested in 
the interplay between the lengthscales of the chemical pattern and the intrinsic lengthscale 
of the contact line. We will therefore consider situations in which the width of the lower 
domains shown in the figure is fixed to ^ ^maxi^e = 90°), while the width of the upper 





23 



domains is varied. In this way, it is expected that the contact hne develops as a finger on the 
lower domains, if the domain is hydrophobic, and spreads sideways to form a sawtooth if the 
domain is hydrophilic. Meanwhile, if the width of the upper domains is sufficiently small, 
the contact line is not expected to deform appreciably. However, for wide enough upper 
domains, a fingering-spreading process of the film, as on the lower domains, is expected. 

Choosing the width of the thinner domains, Wh, as well as the lengths of the hydrophilic 
and hydrophobic domains, Ih and l^, fixes the fraction of the substrate composed of hy- 
drophilic material, fcB = {whlip + w^lh)/[{lh + l4,){wh + w^)]. Given that the widths Wh and 
W(f, are fixed, the limiting cases of = and — oo correspond to patterns of adjacent 
longitudinal stripes, as the ones shown in figures [H](a) and[H](b). In the previous section we 
found that contact line growth is observed in these patterns. Consequently, at intermediate 
values of l^, growth is expected for checkerboard patterns. 

We first consider the case in which the front can destabilize only on the wide domains, 
so we fix Wh to a value that is much smaller than the most unstable wavelength. To gain 
insight on the effect of the checkerboard pattern on the dynamics of the contact line, we 
consider the effect of varying the fraction of hydrophilic domains on the evolution of the 
thin film. We do this by fixing the length of the hydrophilic domains to an arbitrary value, 
which we choose as Ih/xc = 1.3, while the length of hydrophobic domains is varied. The 
width of the upper and lower domains is fixed to Wh/xc = 1.3 and w;^ = 8, respectively. 
We consider six values of l^, namely, l^/xc = 0.6, Itp/xc = 1.0, /^/xc = 1.6, l^p/xc = 2.7, 
l<l)/xc = 4.0 and l^/xc = 8.0. The front is perturbed initially on the wider domains using the 
most unstable wavelength, and on the thinner domains using a wavelength equal to Wh- We 
follow the evolution of the total length of the contact line, L, by measuring its growth rate, 
L = L{wrf),Wh,lti),lh), which is calculated as the difference between the leading and trailing 
edge velocities. The leading edge is located at the middle of the wider domains, while the 
trailing edge is taken as the contact line position at the middle of the thinner domains, as 
shown in figure M 

The evolution of the front is tracked in figures fTOlfa) andfTOlfb) for /^/xc = 0.6, and l^/xc = 
8, which correspond to hydrophilic fractions fcB = 0.25 and fcB = 0.62, respectively. In 
these figures we also show the finger that grows on a homogeneous substrate. For both 
checkerboard patterns, the contact line grows as a finger on a hydrophobic domain, then 
spreads out on a hydrophilic domain, and finally grows again as a finger as it touches the 
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next hydrophobic domain. At a given time, comparing the contact hne for both values 
shows that the leading edge is located at a slightly more advanced position for the smaller 
whereas the trailing edge is located at a slightly less advanced position. This means 
that the leading edge advances faster for small This occurs because the resistance to 
contact line motion is smaller in this case, given that fcB is large. Conversely, at small l^, 
the trailing edge has to sweep increasingly long hydrophobic domains, thus decreasing its 
velocity. Still, whatever the value of Z^, the comparison to the homogeneous substrate shown 
in figure [10] clearly shows that the contact line grows to smaller lengths as a consequence 
of the checkerboard pattern. In figure [I2] we plot the growth rate of the contact line as 
a function of the length of the hydrophobic domains. The limiting cases L{1^ = 0) and 
L{1^ — > oo) correspond to the growth rates of the contact line in the adjacent stripe patterns 
of figures [Ht^b) and[8](a), given that Wh and u;^ are chosen to be equal to the stripe widths 
of those patterns. For intermediate /<^, we observe a monotonous decrease from one value to 
the other. 

Increasing the width of the thin domains, Wh, does not modify the dynamics of the film 
qualitatively, as long as this width is much smaller than Amax- This changes as Wh becomes 
comparable to A^ax] in this case the contact line destabilizes in both domains, giving rise 
to the growth of a different structure. Figure [11] shows the contact line growth on equally 
wide domains, of width equal to A^ax- In this case, there are two fronts that grow with the 
same velocity. Each front grows as a finger when in contact with a hydrophobic domain, 
and spreads on the following hydrophilic domain periodically. At long times, the net growth 
rate of the contact line achieves a constant value. 

Thus, for both symmetric and asymmetric configurations, the checkerboard pattern in- 
troduces the effect of spreading on hydrophilic domains, relaxing the growth of a finger 
momentarily. For asymmetric patterns, this provides a way to tune the growth rate of the 
contact line. In the following section we will propose a kinematical model to account for 
the growth of the contact line in terms of the geometrical parameters of the checkerboard 
pattern. 
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FIG. 12: Contact line growth rate as a function of hydrophobic domain length on checkerboard 
patterns. The growth rate decreases as the length of the hydrophobic domains increases. This effect 
is caused by the adjacent thin hydrophilic domains, whose length increases with l^, thus making 
the trailing edge velocity larger and L smaller. The solid curve corresponds to the theoretical 
prediction, equation ()12p . The straight lines correspond to the limiting growth rates L{1^ = 0) 
(dashed line) and L{1^ — )• oo) (solid line). 

Kinematical model for growth on checkerboard patterns 

To quantify the effect of the checkerboard pattern on the contact line growth, we charac- 
terize the motion of the thin film in terms of the velocities of the leading an trailing edges 
in the wide and thin domains. In general, the instantaneous leading and trailing edge ve- 
locities, vi and Vt, will depend on the configuration of the pattern, i.e., vi = vi{w^, WhJh, l<t>) 
and Vt = Vt{w^, Wh, Ih, l<j>)- We assume that these velocities relax to constant values as soon 
as the contact line comes into contact with a hydrophilic or a hydrophobic domain. For 
simplicity, we take these values as the ones corresponding to the limiting cases of very long 
domains. Therefore, at a hydrophobic domain, the leading and trailing edge velocities are 
approximated as vf = vi{w^, Wh, Ih = 0,1^ ^ oo) and vf = Vt{w^, Wh, Ih oo,l^ = 0), while 
at a hydrophilic domain we have = vi{w^,Wh,lh oo,l^ = 0) and = Vt{w^,Wh,lh = 
0,/^ oo) . 

We denote At; and Atf the times in which the leading and trailing edges sweep a period 
of the underlying pattern, 1^ + It follows that 

At/ = Atf + Atf and At* = Atf + Atf . (10) 

In these expressions, the subscripts / and t, stand for the leading and trailing edges, while 
the superscripts h and stand for hydrophilic and hydrophobic domains respectively. The 
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leading edge will sweep a hydrophilic domain of length Ih with a velocity f f and will then 
move across a hydrophobic domain of length with a velocity vf. On the other hand, 
because of the geometry of the checkerboard pattern, the trailing edge will move over a 
hydrophilic domain of length with a velocity f'*, to later continue over a hydrophobic 
domain of length Ih with a velocity vf. 

Accordingly, equation (fTOjl can be written as 

k±Ih = \ + \ and k±Ill = \ + !L (11) 

Vl vr vr Vt Vt Vt''' 

in terms of the mean values of the leading an trailing edges, vi and Vt, respectively. We can 
finally compute the average growth rate of the contact line a.s L = vi — vt, which within 
approximations considered reads 
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Equation (fT2l) gives an estimate of the growth rate of the contact line in terms of the 
velocities observed in longitudinal striped patterns, the width of the stripes being equal to 
the width of the domains of the checkerboard arrangement. To compare simulation results 
with the kinematical prediction, we measure vf, v^, vf, and v^ from runs that correspond 
to a longitudinal stripe configuration for which Wh and Wrp match with the checkerboard 
pattern, as displayed in figures [8]^a) and[8]^b). Figure [T2l shows a comparison between the 
growth rate measured from simulations and the kinematical prediction as a function of the 
length of the hydrophobic domains. The quantitative agreement shows that the varying 
wettability of alternating patches determines the growth of the contact line essentially by 
fixing the local value of the leading and trailing edge velocities, these values being very well 
approximated by the ones corresponding to long domains. 



V. DISCUSSION AND CONCLUSIONS 



By means of lattice-Boltzmann simulations and kinematical models, we have studied 
the dynamics of driven thin films on a variety of heterogeneous hydrophilic-hydrophobic 
substrates. We have studied the effect of longitudinally stripes and checkerboard patterns. 

We have focused on a scenario where the unstable contact line gives rise to sawtooth and 
finger structures on homogeneous hydrophilic and hydrophobic substrates, respectively, and 
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have examined how the growth of the contact hne is altered by the chemical pattern. To 
this end, we have considered patterns where the typical lengthscale is comparable to the 
most unstable wavelength of the contact line, A^a^. 

On longitudinal patterns, the film follows hydrophilic stripes preferentially, while the 
contact line gives rise to fingering if the width of the hydrophobic stripes is large enough. 
For small enough hydrophobic stripes, the film gives rise to a fingering process caused by 
the longitudinal pattern, where the width of the growing fingers corresponds to the width 
of the hydrophilic stripes. 

On checkerboard patterns, the film undergoes a fingering-spreading process as it moves 
over hydrophobic and hydrophilic domains. Using a kinematical approach, we have shown 
that the net growth rate of the contact line can be tuned by choosing the lengthscale of the 
checkerboard pattern. 

In conclusion, we have performed lattice-Boltzmann simulations of the evolution of a three 
dimensional thin film in contact with chemically patterned substrates We have demonstrated 
that the film can be oriented along a given path, that fingering can be triggered by the 
chemical pattern, and that the growth of the contact line can be controlled by choosing a 
particular configuration for substrate patterning. 
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